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CN ABSTRACT 

Aims. In a follow-up investigation we present Zeeman-Doppler maps of the weak-lined T Tauri star (WTTS) V410 Tau. As a rapid 
rotating star and a typical WTTS the stellar surface of V410 Tau is accessible to surface imaging techniques and allows us to detect 
and reconstruct the major magnetic surface features on this pre-main sequence star. 

Methods. The polarized signals we are measuring are of the order of 10~ 4 to 10~ 3 and are hidden well below the noise level of a single 
observation. A new line profile reconstruction technique based on a singular value decomposition (SVD) allows us to extract the weak 
polarized line profiles (Stokes V) as well as the intensity profiles (Stokes I). One of the key features of the line profile reconstruction is 
that the SVD line profiles are amenable to radiative transfer modeling within our Zeeman-Doppler Imaging code iMap. The code also 
utilizes a new iterative regularization scheme which is independent of any additional surface constraints. To provide more stability a 
vital part of our inversion strategy is to invert both Stokes I and Stokes V profiles to simultaneously reconstruct the temperature and 
magnetic field surface distribution of V410 Tau. A new image-shear analysis is also implemented to allow the search for image and 
C"| line profile distortions induced by a differential rotation of the star. 

Results. The magnetic field structure we obtain for V410 Tau shows a good spatial correlation with the surface temperature and is 
dominated by a strong field within the cool polar spot. The Zeeman-Doppler maps exhibit a large-scale organization of both polarities 
around the polar cap in the form of a twisted bipolar structure. The magnetic field reaches a value of almost 2 kG within the polar 
region but smaller fields are also present down to lower latitudes. The pronounced non-axisymmetric field structure and the non- 
detection of a differential rotation for V410 Tau supports the idea of an underlying a 2 -type dynamo, which is predicted for weak-lined 
T Tauri stars. 
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1 . Introduction fortunately, not only true for WTTS stars but also for any other 

r^sl type of star (Schrijver & Zwaan 2000). Even for the Sun there 

- ^. T Tauri stars are very young stars on their way towards the zero- h nQ dear answer of how strQng sunspQts and actiye regions 

-H age main sequence (ZAMS). With the ongoing contraction they are stm connected d down int0 the tachocline reg ion (Rempel 

prepare the conditions in their interiors to start the first nuclear & Schussler 2009; R j 201 1} where the internal d is 

process. In classical T Tauri stars, the outside is still dominated believed tQ operate (Schrijver & xitle 1999; Sc hussler 2005). 

by the accretion of matter from the disk while in the weak-lined We ^ m in the hase of characterization an d it may still take 

T Tauri stars (WTTS), or naked T Tauri stars, the accretion SQme dme before wg C£m connect observations with theory (and 

disk is basically gone and the star is soon to arrive on the ZAMS. yice versa ) j n a conclusive way 
, Little is known about the magnetic field generating processes TT ' .. . 

b during both of these evolutionary phases (see, e.g. the review by However the characterization of surface magnetic fields with 
3 Gregory & Donati 201 1). How do WTTS stars produce, main- llke Do PP f ^ging (DI (e.g. Rice 2002) which 

tain, and reorganize their internal fields at this rapidly changing tracks the magnetically induced cool spots on the surface of 

evolutionary stage? An even more relevant question from a gen- fare, or more direct y with Zeeman-Doppler Imaging (ZDI) 

, • . e ■ ■ , j c , , - t . , (see, e.g. Donati 2001), has provided us in the last two decades 

eral point of view is how do these fields manifest themselves v . 6 . . v . 

on the surface of the star and whether their appearance, i.e. the wlth a weal * of new lnfo ™ n about surface s P° ts and 

topology, is directly linked to the generating process in the inte- magnet ^ n ^ lds °" many d \ T type ' °L S , tarS ' See StraSS " 

rior? Or do convection and other surface flows alter the surface meier ^° 09 ) and Donatl & Landstreet ^ 2009 ) for comprehen- 
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appearance in such a way that they will prevent us from making 

any direct conclusions about the the internal dynamo processes? In this work , we concentrate on DI and ZDI of V410 Tau, a 

We are still far from answering all these questions which is, un- ra P idl y rotating WTTS that has fully stripped away its surround- 

ing disk. Many WTTS like V410 Tau show strong signatures 
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activity rather than to a left-over from external accretion pro- 
cesses. Indirect signs of magnetic activity are known for a long 
time from photometric light curve variations (see, e.g. Rydgren 
& Vrba 1983; Bouvier & Bertout 1989). More direct evidence 
of photospheric magnetic activity could be obtained by Doppler 
images, in particular for V410 Tau (Joncour et al. 1994; Strass- 
meier, Welty & Rice 1994; Hatzes 1995; Rice & Strassmeier 
1996; Rice, Strassmeier & Kopf 2011). Skelly et al. (2010) 
presented the first Zeeman-Doppler map of V410 Tau. Surpris- 
ingly, the reconstructed large-scale field exhibit a rather complex 
structure dispersed over different latitudes with strongly inclined 
fields that are only weakly correlated with the spots shown in 
their brightness maps. 

Signs of chromospheric and coronal activity on V410 Tau 
are also known for many years and strong ultraviolet and opti- 
cal emission lines and X-rays emission were observed (Hatzes 
1995; Stelzer & Neuhiiuser 2001; Stelzer et al. 2003). An in- 
triguing fact here is that the modulation of the chromospheric 
and coronal activity tracers seem not to be strongly correlated 
with the photospheric spot locations (Stelzer et al. 2003). 

To shed more light on the magnetic field distribution on stel- 
lar surfaces, in particular for V410 Tau, we conducted a simulta- 
neous investigation of rotationally modulated spectral line pro- 
files (Stokes I) and spectropolarimetric line profiles (Stokes V) 
in terms of a rigorous Doppler imaging and Zeeman-Doppler 
imaging analysis. In this respect, the present paper is a contin- 
uation of the work by Rice, Strassmeier & Kopf (201 1), which 
we hereby call Paper I. The present paper is organized as fol- 
lows. In Sect. 2, we briefly describe the acquisition of the data. 
In Sect. 3, we introduce our improved polarized-line profile ex- 
traction and reconstruction technique based on a singular-value 
decomposition (SVD). Sect. 4 is devoted to a recap of the ba- 
sic outline and assumptions in our ZDI code iMap and its recent 
improvements. In Sect. 5, we present the results of our simulta- 
neous DI and ZDI approach as well as the results of our search 
for a differential rotation on V410 Tau. Sect. 6 summarizes our 
findings and conclusions. 

2. Observations 

High-resolution spectropolarimetric observations were obtained 
with the ESPADONS echelle spectrograph and polarimeter (Do- 
nati 2003) in Stokes I and Stokes V at the 3.6m Canada-France- 
Hawaii telescope (CFHT) on Mauna Kea, Hawaii. Data were 
obtained in queue mode during one observing block over four 
nights from 2008 October 15 to 19, and one over 13 nights dis- 
tributed between the dates of 2008, December 5 to 2009, Jan- 
uary 14. A spectral resolution of 60,000 with a wavelength cov- 
erage of 390-900 nm was obtained. All integrations on V410 
Tau (V ~ ll m .0) were set to an exposure time of 4 x 600 s. 
This allowed for a total of 17 spectra with an average signal-to- 
noise (S/N) ratio of around 160:1 per single exposure per pixel. 
All spectra were reduced and extracted using the Libre-ESPRIT 
package provided by CFHT and executed automatically. Details 
of this reduction procedure are given in (Donati et al. 1997). 

3. SVD based signal extraction and reconstruction 

3.1. Method 

For most active cool late-type stars the typical signal amplitude 
of a single polarized line profile is of the order of a few times 
10 4 relative to the normalized continuum. These values are far 
below the signal-to-noise ratio of a typical spectropolarimetric 



observation, which hardly exceeds a few times 10 3 . So the indi- 
vidual polarimetric line profiles are buried well below the noise 
level by typically one order of magnitude. This requires tech- 
niques that are able to recover and extract the line profile sig- 
nals from the measurement noise prior to the actual ZDI analy- 
sis. Most of the techniques for circularized, i.e. Stokes V, line 
profiles take advantage of the fact that the polarized Stokes V 
patterns for almost all Zeeman sensitive atomic spectral lines 
are similar except for a non-linear scaling (in wavelength and 
amplitude) caused by different atomic line parameters, broad- 
ening mechanisms, and magnetic sensitivities. For rapidly ro- 
tating stars, when the rotational broadening dominates as well 
as under the weak magnetic field approximation (Stenflo 1994), 
line-profile shapes can be considered as equal and amplitudes are 
only linearly dependent upon their individual Lande factors and 
line strengths. One popular and very successful extraction tech- 
nique which follows this line of arguments is the so called least- 
squares deconvolution (LSD, Donati et al. 1997; Kochukhov et 
al. 2010). Another method is the principal component analysis 
(PCA, Carroll et al. 2007; Martinez Gonzalez et al. 2008; Car- 
roll et al. 2009; Paletou 2012) or the simple but very effective 
coherent addition of line profiles in the velocity or logarithmic 
wavelength domain (Semel et al. 2009; Ramirez et al. 2010). 

For our application to V410 Tau, we present here an im- 
proved version of the PCA technique that is more robust in the 
case of weak polarimetric signals and which can be described in 
terms of an eigenvalue decomposition (EVD) of the observation 
covariance matrix. An EVD of the signal covariance matrix (or 
a SVD of the observation matrix) allows us to obtain an orthog- 
onal decomposition of the signal or line profiles in terms of the 
covariance eigenprofiles. The basic idea is similar as in the pre- 
vious PCA-based approach where the similarity of the individual 
Stokes V profiles allows one to describe the most coherent and 
systematic features present in all spectral line profiles as a pro- 
jection onto a small number of eigenprofiles. In parallel, inco- 
herent features like noise, and line blends etc. will be dispersed 
along many dimensions in the transformed eigenspace. In other 
words the similarity and coherence of all individual Stokes V 
patterns allows one to capture the signal pattern by a low-rank 
representation of the observation covariance matrix. The respec- 
tive eigenvalues of the observation covariance matrix, together 
with an estimate of the noise level in the data, then provide a 
means of separating the low-dimensional signal subspace from 
the noise subspace. Once the reduced rank of the signal sub- 
space is determined a projection of the observed line profiles 
onto the signal subspace and a subsequent subspace-averaging 
is performed to eventually provide the signal boosting effect for 
the polarimetric line profile. What is described in the following 
for Stokes V profiles also holds for Stokes I but for the sake of 
brevity we will use a notation that makes use of Stokes V only. 

The problem setup is as follows: An individual, observed, 
line profile V(/l) is the result of the true signal vector V(A) and 
a noise vector N, which we consider as additive and indepen- 
dent of wavelength and time and of zero mean. Although our 
idealize assumption are not met (and we will discuss the conse- 
quence later in this section), let us proceed with this assumption 
to describe the basic mechanism of the subspace method. 

The observed line profile can be written in a discrete form 

as, 

% = Vij + N , (1) 

where i and j denote the wavelength and the spectral line index, 
respectively. In vector form, we may write Eq. 1 as 

V,- = V,- + N . (2) 
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We further assume that the net circular polarization of each indi- 
vidual Stokes V line profile is zero such that their integral over 
the wavelength range is also zero. This implies that there are 
no gradients in magnetic field and velocity present in the atmo- 
sphere (Carroll & Kopf 2007) which is also the usual starting 
point for ZDI to keep the problem tractable. 

Arranging now the individual, observed, Stokes profiles in 
a m x n observation matrix V, where m denotes the number of 
wavelength pixels and n the number of observed individual spec- 
tral lines, we can calculate the non-normalized covariance matrix 
of the observations as 



Cy = W 



(3) 



where the superscript T denotes the transpose. The covariance 
matrix Cy is a ra x ra real symmetric matrix and has m linear 
independent eigenvectors, which can be used to decompose Cy 
by a spectral decomposition, i.e. an eigenvalue decomposition, 
such as 



QAQ 



(4) 



with A being a diagonal matrix containing the eigenvalues i„ 
and Qanraxm orthogonal matrix containing in each column 
the eigenvectors q,. By assumption, signal and noise are uncor- 
rected and we may therefore write the covariance matrix C v as 
the sum of the two matrices Cs and Cn, i.e. the clean signal and 
the noise covariance matrix: 



Cs + Cn . 



(5) 



Based on our assumption that the clean signals span only a lower 
dimensional subspace of the observed space (i.e. signals plus 
noise), which is equivalent of assuming that the signal matrix 
is rank deficient, we can use the orthogonal matrix Q of the 
observed signal space for the eigenvalue decomposition of the 
clean signal covariance matrix Cs and the noise covariance ma- 
trix Cn, to write 



and 



Q A s Q 1 



Q (cr 2 !) Q 1 



(6) 



(7) 



which directly follows from Eq. (5). Here, a 2 is the eigenvalue 
of the noise and corresponds to the noise variance along a direc- 
tion in subspace, the matrix I is the identity matrix and A s is the 
diagonal matrix of the eigenvalues A of the clean signal matrix. 
Note, because the noise is assumed to be isotropic and feature- 
less it will be homogeneously distributed across the eigenspace 
with all eigenvalues of similar magnitude. Note also that by defi- 
nition of the covariance matrix each eigenvalue provide the total 
energy (variance) of the observation along its corresponding di- 
mension (eigenprofile). The noise eigenvalue is related to the 
estimated observation error cr ofcs by 



2 2 
cr — m cr 



obs 



(8) 



Using the definition of Eqs. (6) and (7) the eigenvalue decompo- 
sition of the covariance matrix C v can now be written as 



CV = Q(As +cr 2 I)Q r . 



(9) 



If the clean signals are confined to an s-dimensional subspace 
with s < m, then the eigenvalues of A can be expressed as 



At + cr 2 
cr 2 



i = 1...S 

i — s + 1...W . 



(10) 



Partitioning Q in an orthogonal m x s matrix Qs for the signal 
subspace and a m x (m - s) matrix Q N for the noise subspace 
allows us to write for the observational covariance matrix 



[QsQn] 



A s 




+ cr 



I„ 







[QsQnY 



(11) 



From this we can define a simple but effective criterion for the 
separation of the signal subspace. All those eigenvalues A t , i.e. 
the combination of clean signal and noise eigenvalues, which are 
larger than the noise eigenvalues alone carry more energy than 
a pure noise component and may therefore exhibit significant 
structural information. The condition for the separation of the 
signal space then reads 

A t > cr 2 . (12) 

It provides an estimate of the number s of the signal space, i.e. 
the reduced rank. The above described eigenvalue decomposi- 
tion can efficiently calculated by utilizing a singular value de- 
composition (SVD) of the original observation matrix V. 

Let us now briefly discuss the restriction and consequences 
of our idealized assumption, i.e. that we deal with homogeneous 
and stationary noise. The real noise that we are dealing with 
in our spectral data is certainly not stationary and may have a 
dependence on wavelength. The noise can therefore introduce a 
systematic feature and the noise matrix will have off-diagonal el- 
ements. If the correlation of this feature is strong enough among 
the contributing lines it would give rise to a structure in signal 
space that is different from pure random noise. In that case this 
systematic will therefore be detected as a significant signal fea- 
ture in our framework as it follows our basic notion of being 
systematic and correlated among all contributing lines. So with- 
out giving up the derivation of the rank estimation we need to 
change the meaning of the signal space. The signal space is not 
exclusively given by the true underlying Stokes signal but rather 
includes all systematic effects that possess a certain correlation 
within the observation matrix. 

In principle, we could now retrieve each individual Stokes V 
profile by projecting the observed line profiles V,- onto the re- 
stricted set of s signal eigenprofiles from Q s . This could be 
done for any line profile of interest. Unfortunately, this approach 
alone is limited to relative noise levels above unity. The situa- 
tion changes when the original profiles are noise dominated such 
that the true signal patterns are completely hidden by noise and 
any projection operations is increasingly affected by the noise 
contamination of the original signal. For this case we use a two- 
stage strategy where the reduced rank estimation of the covari- 
ance matrix is used to perform a coherent addition of the pro- 
jection coefficients in signal subspace. Consequently our second 
stage of the signal extraction consists of a projection of each in- 
dividual line profile onto the eigenprofiles Qs in order to obtain a 
low-dimensional representation of each observation vector. This 
projection and the respective matrix of projection coefficients Vp 
can be written as 



Vp = Qs V 



(13) 



We can use the projection matrix to obtain the subspace aver- 
aged projection vector Vp by multiplying an n x 1 sum vector 
Z to the projection matrix Vp. Each entry of the vector Z may 
hold a specific weighting for each individual spectral line or may 
be given by a simple equal weighting scheme (e.g. 1/n for the 
sample mean). 



Vp = VpZ. 



(14) 
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We can then use the averaged projection vector V P together with 
the signal space eigenprofiles to finally obtain the signal recon- 
struction Vr, 



(15) 



Now we have a reconstructed (weighted) mean line profile that is 
comparable to the mean Zeeman feature obtained by a conven- 
tional LSD approach. However, the reconstruction here is per- 
formed over the estimated signal subspace, which carries the ma- 
jority of the significant information of the profile pattern. This 
is in contrast to the LSD method which operates in the origi- 
nal data domain where the observed signals are still a mixture 
of coherent and non-coherent features. Much of the noise and 
incoherent structures are already reduced in the first process- 
ing step while the main signal-boosting is provided by Eq. (14). 
The second crucial difference lies in the fact that we use this 
subspace-averaged line profile not as a proxy for the general line 
shape in the ZDI inversion but rather as the real approximate 
mean line profile that can be synthetically calculated by using 
the same lines that have been used for the observation matrix V. 
Even though it is a formidable computational task to calculate 
the mean line profile from several hundreds or even thousands 
of individual spectral lines in each iteration cycle of the inver- 
sion process it has been shown that artificial neural networks 
(Bishop 1995) provide a fast and accurate approximation of the 
polarized radiative transfer equation to calculate thousands of 
synthetic Stokes profiles with a minimum amount of computa- 
tional time (Carroll, Kopf & Strassmeier 2008). 

The net gain of the SVD method can be readily estimated 
by using our starting assumption of incoherent white noise that 
is homogeneously distributed in eigenspace. The total energy 
(variance) of a line profile is given by the following expectation 



o-l = E{\ T \} = £ \ . 



(16) 



We know that the reconstructed signal only uses eigenprofiles of 
the signal subspace so that we can write according to Eq. (11) 
their total energy as, 



(17) 



i=l 



The relative signal-to-noise ratio (SNR) of the reconstructed sig- 
nal can then be expressed as 



SNR S 



(18) 



On the other hand using again Eq. (10) as well as Eq. (16) the 
energy of the original signal which results from the signal sub- 
space eigenvalues plus the noise contribution from the entire tri- 
dimensional domain is given by, 



<=i 



and the relative signal-to-noise ratio in this case is 



SNR m = 



r2 • 



(19) 



(20) 



This finally allows us to write the gain factor g for the first stage 
of the SVD reconstruction, 



9 = 



SNR S 
SNR„ 



m 
s 



(21) 



We see that due to the projection of the observations onto the 
signal subspace we already obtain a signal improvement of the 
order of m/s, i.e. the noise level scales with -\fmjs. The sample 
size also plays an import role in the reconstruction process and 
is relevant for the noise level within the signal eigenprofiles as 
well as for the second stage of the reconstruction process. Using 
Eq. (15) together with the mean weighting scheme in the second 
stage and moreover realizing that each subspace dimension con- 
tributes a total noise variance which is equal to the eigenvalue 
cr 2 , the variance of the average of all subspace dimensions then 
scales with cr 2 /n. 

However, let us again emphasize that our entire derivation 
above is based on an idealized noise behavior (i.e. wavelength 
independent white noise) which is might be not fulfilled in real 
high-resolution spectroscopic observations. For a more reliable 
rank estimation as well as an estimation of the noise level for 
the reconstructed Stokes profiles, we will rely on a randomized 
bootstrap procedure (see next section). 



3.2. Application toV410 Tau 

To prepare the observation matrix V, we first have to select a 
suitable number of Zeeman-sensitive spectral lines. This is done 
by calculating a synthetic spectral atlas to judge which spectral 
line profiles are strong enough and provide enough information. 

By utilizing the forward module of iMap (Carroll, Kopf & 
Strassmeier 2008) we synthesize Stokes I and Stokes V spec- 
tra for the observed wavelength range of our data (480-740nm). 
The individual line parameters are taken from the VALD line list 
(Piskunov et al. 1995; Kupka et al. 1999). The astrophysical pa- 
rameters for the synthetic calculations are those summarized in 
Sect. 4.4, except for the rotational velocity which is set to zero. 
A uniform temperature star with a homogeneous and purely ra- 
dial surface magnetic field of 500 G is taken for the model spec- 
tra. The model atmospheres used throughout the analysis are 
Kurucz/Atlas-9 models (Castelli & Kurucz 2004). 

The resulting synthetic spectral atlas (line list) in Stokes I 
and Stokes V is then used to select those spectral lines that fulfill 
the following criteria: Stokes I line depth > 0.6 and Lande factor 
> 1.5. 

This procedure allowed us to select a list of 929 spectral 
lines. The majority of these line are neutral Fe, Ca, Ti, Na, and C 
lines. For Stokes I, we only select a short list of the strongest 56 
iron lines because the intensity profiles are two to three orders of 
magnitude larger than the polarized spectra. In this case there is 
no need for a full reconstruction but instead a denoising suffices, 
which requires a significantly smaller number of spectral lines. 
In a first step, the observed spectra are transformed into velocity 
domain with a step size of 1.5 kms -1 , which corresponds ap- 
proximately to the lowest wavelength sampling in the observed 
spectra. Then the line list is used to extract the corresponding 
Stokes I and Stokes V profiles within a range of ±150 kms -1 
around their respective line centers. Using the extracted line pro- 
files, we then build the Stokes V observation matrix V and the 
Stokes I observation matrix I as described in Sect. 3. 

In total, 12 rotational phases were obtained within a time 
span of four weeks during December 2008 and January 2009. 
For each of these 12 rotation phases, indexed with p, we create 
an observation matrix V ; , and l p with 201 variables (intensities 
along the velocity domain) and a sample size of 929 for and 
56 for lp. As described in Sect. 3, a singular-value decompo- 
sition provides us with the eigenvalues and eigenprofiles of the 
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Fig. 1. Stokes V eigenvalue spectrum of a V410-Tau observation (solid 
line) and the mean eigenvalue spectrum obtained from a bootstrap ran- 
domization process (dashed line) for the first 100 dimensions of rotation 
phase 0.014. The magnitude of the eigenvalues of the corresponding ob- 
servation covariance matrix and the mean eigenvalues of all randomized 
matrices are plotted over the corresponding subspace dimension index. 
The crossing point at dimension 27 marks the dimension from which on 
the individual eigenprofiles contain no more significant information. 

0.0006 f i 




Subspace dimension 



Fig. 2. Normalized Stokes I eigenvalue spectrum of the observation 
(solid line) and the noise variance (dashed line) for the first 20 subspace 
dimensions of rotation phase 0.014. For the purpose of a better illustra- 
tion, the first three eigenvalues are not shown since they are significantly 
larger. The crossing point determines the signal subspace dimension to 
be 5. 



observation covariance matrix from which we derive for each 
phase a reconstructed average SVD profile. 

Because the original raw Stokes V spectra show no apparent 
sign of a detectable line profile polarization, i.e the polarimet- 
ric signals are much weaker than the noise level, we adopt as an 
initial hypothesis that the observed Stokes V spectra contain no 
signal information except pure noise (not necessarily white). If 
this hypothesis is true then a random permutation of the intensity 
values of each line across the velocity domain should not change 
the outcome of our SVD analysis. In other words, if there are 
no systematic and coherent features present among the spectral 
line profiles in the observation matrix, a random scrambling of 
their individual wavelength or velocity binned intensity values 
would give rise to the same series of eigenvalues. This series 
of ordered eigenvalues (ordered according their magnitude) is 
called the eigenvalue spectrum. We therefore build a second ob- 



servation matrix for Stokes V where each extracted line profile 
is subject to randomization across its velocity domain. If there 
are some sort of systematic and coherent structures present in all 
observations the randomization process will eliminate them such 
that the randomized covariance matrix provides essentially a flat 
spectrum of eigenvalues, i.e. all eigenvalues A* have the same 
value A\ - A* 2 - ... - A* m . Hence, in the case where no significant 
systematic structures are present in the original observation ma- 
trix the eigenvalue spectrum of the observation and randomized 
covariance matrix will result in a flat eigenvalue spectrum within 
some error margin. On the other hand if there are correlated and 
systematic features present we expect a distinct difference in the 
eigenvalue spectrum. 

However, from random matrix theory it is known that a pure 
random matrix does not necessarily provide a flat eigenvalue 
spectrum and depend on the number of variables (wavelength) 
and measurements (spectral lines) (see e.g. Johnstone 2001). 
For this reason, we perform a Monte-Carlo simulation to provide 
a large sample (1,000) of randomized covariance matrices from 
which we can deduce the mean eigenvalue spectrum and robust 
error margins for each eigenvalue. The errors of the eigenval- 
ues are estimated by a bootstrap procedure (Efron & Tibshirani 
1994) where we used again a sample size of 1,000. After calcu- 
lating the SVD of all the randomized matrices, we use a parallel 
analysis (Jolliffe 2002) to compare the significance of the eigen- 
values of the observation covariance matrix relative to the mean 
eigenvalues (averaged over the respective subspace dimension) 
of all randomized covariance matrices. 

Note, that all matrices, the original observation matrix and 
each of the randomized matrices, carry the same amount of total 
energy (variance), i.e. the total sum over all of the eigenvalues 
are equal. After sorting the eigenprofiles for all covariance ma- 
trices according to the magnitude of their respective eigenvalues, 
we can finally identify the eigenprofiles of the observation co- 
variance matrix that contain significant profile information. This 
means we determine all those eigenprofiles whose eigenvalues 
A, are larger then the corresponding mean eigenvalues (A*) of 
the randomized matrices. This provides a noise threshold value 
according to Eq. (12) such that A t > {A*} = a 2 which in turn 
allows us to determine the effective rank of the signal subspace. 

Fig. 1 shows the eigenvalue spectrum of the observation ma- 
trix and the mean eigenvalues from all the randomized matrices 
for one observational phase (0.014). The magnitude of the lead- 
ing eigenvalues (solid line) obtained from the observation matrix 
provides a first evidence that the observation matrix contains sig- 
nificant information different from pure noise or random fluctua- 
tions which is indicated by the randomized eigenvalue spectrum 
(dashed line). The crossing point at which the eigenvalues from 
the observation covariance matrix drop below the eigenvalues of 
the randomized matrices marks the critical dimension in the sub- 
space. Up to eigenvalue 27 the respective eigenprofiles contain 
significant information above the noise level. The estimated re- 
duced rank is thus 27. This estimate allows us to proceed and 
use Eq. (15) to calculate the rank-reduced (subspace) averaged 
Stokes V profile. This procedure is then performed for all rota- 
tion phases to reconstruct the entire set of Stokes V profiles used 
in the following DI/ZDI inversion. Note that we omitted the er- 
ror bars for both eigenvalue spectra in Fig. 1 because the standard 
deviation for all values is lower then 2% , which hardly exceeds 
the thickness of the plotted lines. For the Stokes I profiles, we 
have a significant larger relative signal-to-noise ratio (it is in fact 
the S/N estimate from 2). Therefore, we do not start from a neg- 
ative hypothesis that no signal or coherent information is present 
in the data set but instead directly apply Eq. (12) to determine 



Article number, page 5 of 18 



0-0020 
0.0015 
0.0010 
0-0005 
0.0000 



Phase A 0-014000 
Std. Error : 0-000102. 




-150-100 -50 50 100 150 
Velocity [km/s] 




-150-100 -50 50 100 150 
Velocity [km/s] 



I -50 50 100 150 
Velocity [km/s] 



Phase : 0.251000- 
Std. Error : 0-000104 - 




-150-100 -50 50 100 150 
Velocity [km/s] 




-150-100 -50 50 100 150 
Velocity [km/s] 



hase : 0-341600 
td. Error : 0.000102" 




-150-100 -50 50 100 150 
Velocity [km/s] 



I -50 50 100 150 
Velocity [km/s] 



'Jfcse : 0.41 1300 
SlA Error : 0.000124 




-150-100 -50 50 100 150 
Velocity [km/s] 




-150-100 -50 50 100 150 
Velocity [km/s] 




-150-100 -50 50 100 150 
Velocity [km/s] 



) -50 50 100 150 
Velocity [km/s] 



Fflhase : 0.906200 
ltd. Error : 0.000128 




-150-100 -50 50 100 150 
Velocity [km/s] 



Fig. 3. Reconstructed SVD Stokes V profiles for each rotational phase of V410 Tau. Each profile is reconstructed from an observation matrix 
with a sample size of 929 spectral lines . The averaged standard error is given together with the corresponding rotational phase. 
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Fig. 4. Reconstructed SVD Stokes I profiles for each rotational phase. For all reconstructed profiles an observation matrix with 56 spectral lines 
is used. The averaged standard error is given together with the corresponding rotational phase. 



the signal subspace. Fig. 2 shows the normalized eigenspectrum (solid line) and the error variance o- a t s as determined from the 
for the first 20 eigenvalues of the observation covariance matrix observation (dashed line). It is interesting to see that almost the 
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entire energy (variance) of the observation covariance matrix can 
be described by a very small number of eigenprofiles. The cross- 
ing point indicates that the subspace dimension is five. The first 
five eigenprofiles contain already 99.98% of the total energy of 
all line profiles within the observation matrix. In the following 
step, we again use Eq. (15) to calculate the subspace averaged 
Stokes I profiles. 

Finally, Fig. 3 and Fig. 4 show the reconstructed set of 
Stokes V and Stokes I profiles for all rotational phases. The 
standard errors are provided again by a Bootstrap method with 
a re-sampling number of 1,000. For each velocity bin, we thus 
obtained an estimate of the standard error which is condensed in 
Fig. 3 and Fig. 4 to a mean standard error value averaged over 
the velocity domain. Its typical value is 1 x 10~ 4 for Stokes V 
and 1 x 10~ 3 for Stokes I which corresponds to a signal-to-noise 
level of 1 0,000 or 1 ,000 respectively 



4. Stellar surface reconstruction with iMap 

4.1. A brief description of iMap 

In this study, we use the iMap code (Carroll et al. 2007, 2009) 
to simultaneously reproduce the temperature and magnetic vec- 
tor field distribution from a sequence of observed Stokes I and 
V profiles. The forward modeling in iMap is based on polarized 
radiative transfer to allow for the best possible accuracy in line 
profile modeling (Carroll et al. 2009). All line parameters nec- 
essary for the synthetic calculation are taken from the Vienna 
atomic line list (VALD) (Kupka et al. 1999). The underlying 
model atmospheres are based on Kurucz ATLAS-9 atmospheres 
(Castelli & Kurucz 2004). The actual calculation and integra- 
tion of the polarized radiative transfer equation in LTE can be 
performed either by a numerical integration or by a number of 
trained artificial neural networks (ANN) (Carroll, Kopf & Strass- 
meier 2008). 

The iMap code is equipped with a new inversion module. 
While the former versions relied on a conjugate gradient method 
(Press et al. 1992) with a local entropy regularization (Car- 
roll et al. 2007), the current version of iMap uses an itera- 
tively regularized Landweber method (Engl, Hanke & Neubauer 
1996). Many different regularization methods exist for linear 
and non-linear inverse problems but interestingly only two of 
the most restrictive ones received attention in the field of ZDI 
and DI namely the Tikhonov regularization (see e.g. Piskunov & 
Kochukhov 2002) and the maximum entropy method (e.g. Vogt 
et al. 1987). There is also a long list of comparisons between 
these two methods and their benefits and respective shortcom- 
ings (see, e.g. Rice 2002). The three major issues that constitute 
an ill-posed, inverse problem are : (i) a solution may not exist, 

(ii) the problem may not be unique, (iii) the solution does not de- 
pend continuously on the data. We will not go into problem (i) 
which would require a discussion about the adequacy of the un- 
derlying ZDI/DI model assumptions (see Carroll et al. 2009, for 
a discussion), but instead we assume (as usual) that our model 
admits a solution. Then, we are left with the uniqueness prob- 
lem (ii) and the numerical stability problem (iii) which both de- 
pend also on the noise level of the data. Problems (ii) as well 
as (iii) are closely related and can effectively be addressed by 
regularization methods that put additional constraints on the so- 
lution. Please also see in this context the paper of (Kochukhov 
& Piskunov 2002) who describes the degeneracy caused by us- 
ing just the information of the Stokes I and Stokes V profiles. In 
our iterative approach, implemented in iMap, problem (ii) and 

(iii) will be addressed by an iteration process where the step size 



and an appropriate stopping rule provide the regularization of 
the inverse problem (Engl, Hanke & Neubauer 1996). 

Iterative regularization for inverse problems has been the 
subject of various theoretical investigations over the recent years 
(Hanke 1997; Engl & Kugler 2004; Egger & Neubauer 2005; 
Kaltenbacher et al. 2008). The Landweber iteration, which is 
used here, rests on the idea of a simple fixed-point iteration, de- 
rived from minimizing the sum of the squared errors. Our new 
inversion routine follows exactly this line and can be described 
as follows: Written in a concise vector notation the problem set- 
ting is 



min ^||I(x) 

x 2 



■on 2 



(22) 



where ||| is the L2 norm and I is the synthetic model profile 
over all spectral lines, wavelengths or velocities, and rotational 
phases, O is the corresponding observation. The vector x con- 
tains all our free parameters of the model, i.e. the temperature 
and the magnetic field vector for each surface element. The it- 
eration now proceeds along the negative gradient direction and 
updates the current estimate of the solution vector, x k , in the fol- 
lowing manner 



Xjt+i = x k + w k T(\ k ) (O - I(Xjfc)) . 



(23) 



Here, I' is the gradient vector with respect to all surface element 
values and w k is the weight factor that can adaptively accelerate 
the iteration process. In the conventional Landweber iteration 
process, w k is set to unity. To accelerate the procedure we use a 
variant of the steepest descent (Kaltenbacher et al. 2008) and set 
w k to 



w k = 



INI 2 



P'(x*)w*|l 2 



(24) 



where u k = I' * (O - l(x k )). 

The semi-convergence (Hanke, Neubauer & Scherzer 1995) 
of the method requires a stopping rule before it enters into the 
noise level of the data to regularize the procedure and to avoid 
overfitting. One common, and well studied criterion for the stop- 
ping condition is the Morozov discrepancy principle which can 
be written for the iterative approach as 



Hlftr) - 0|| < t5< Hlfct) - 0|| ; < k < k* , 



(25) 



where 6 is an upper bound for the data error (i.e. noise), r a 
positive number and k* the maximum iteration number. In terms 
of stability and convergence it can be shown that r has to satisfy 
t > 1 (Engl, Hanke & Neubauer 1996). Formally the noise esti- 
mate for the ZDI/DI can be derived from the observations. Given 
a noise contribution <r, for each velocity i we can combine the in- 
dividual errors to form a vector <j such as, cr T = (cri , <x 2 , cr„) 
with n the number of velocity points. The error can then be ex- 
pressed as 6 - ||cr||. If the noise estimate is homogeneous (i.e. 
equal for all velocities) or does not vary much over the velocity 
domain and time we may use the maximum of <x to write for 6 
the relation 



5 = max(cr) yfn . 



(26) 



If we assume that Eq. (22) follows a x 1 distribution we may use 
the number of degrees-of-freedom of the problem to determine 
the error resulting from the limited degree of freedom of the 
model to write 



5 = max(cr) yjn - p. , 



(27) 
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where p is the number of parameters in the model. If the inver- 
sion reduces the error function Eq. (22) down to the threshold 6 
we would ensure that the reduced x 2 is close to one. But let us 
emphasize here that DI as well as ZDI (given we have correctly 
modeled the problem) are non-linear problems (see Carroll et al. 
2009) and moreover the problem is generally ill-posed, which 
makes it by no means a simple task to adequately determine the 
real degree of freedom of the problem. The reduced^ 2 as a mea- 
sure of the goodness-of-fit may therefore only of limited use. 

If the data are pre-processed like in this work the noise esti- 
mate for the SVD reconstructed data can be directly obtained 
from the error estimates provided by the bootstrap procedure 
given in the last section. 

The simultaneous DI and ZDI procedure for each surface el- 
ement, is implemented in the following way. Instead of using 
one long vector x that incorporates the entire information of the 
surface temperature and the magnetic field, we proceed with an 
alternate minimization (Byrne 2004) of Eq. (22). Alternate min- 
imization means here that after each iteration for the temperature 
x^ +1 a parallel inversion of the magnetic field vector x B is per- 
formed, which uses the information of the k + 1 iteration of the 
temperature minimization. The alternate iteration at step k + 1 
for x^ +1 and x^ +1 can then be written as 

x T M = x[ + ^I'(x[,xf)(0-I(x[,xf)) , 

xf+i = 4 + . 4) (o - i(x[ +1 , xf )) . 



(28) 



provides a reconstructed composite line profile. The SVD re- 
construction process extracts a single (weighted) mean SVD 
line profile that comprises the information of all the contribut- 
ing lines of the observation matrix and therefore the error func- 
tion can no longer be applied to the sum of squared differences 
of the individual lines but needs to be applied on the squared 
differences between the reconstructed SVD profile and a proper 
synthetic equivalent. The strategy here is in fact to compute a 
synthetic mean profile for all contributing spectral lines used in 
the SVD reconstruction process by accounting for all radiative 
transfer effects and line modeling characteristics (e.g. model at- 
mospheres, spectral line parameters and line blends). 

We briefly show that the result of this strategy provides an 
error which is on average at least as low as in the conventional 
case, where the sum of differences from each contributing line 
is computed according Eq. (22). Let us begin with the conven- 
tional case of Eq. (29). For brevity we incorporate the rotational 
phases and the wavelengths or velocities into the vector notation 
such that \k and Ok represent the Stokes profile of one particular 
spectral line k taken over the entire velocity range and over all 
available rotational phases. If we correctly model each individ- 
ual spectral line that contributes to the averaging process, then 
we can express the difference between the observed Stokes pro- 
file O/i and the synthetic profile I* as being exclusively due to 
the inherent noise et such that I* corresponds to the observation 
Ok within the limits of e^; 



I* 



Note that both x^ +1 and x* +1 depend on each other due to radia 
tive transfer effects. 

We want to note here that after a careful analysis of the con 
vergence and stability behavior of the inversion, the strict simul- 
taneity of the temperature and magnetic inversion is not neces- Ek - & [(It - O^) 2 ] = & |e 2 ] 
sarily required whereas the order of the inversion is of impor- 
tance. The dependence of the DI on the magnetic inversion is not 
strong, this is mainly due to the weakness of the effect that the 
magnetic field introduces via the Zeeman effect on the intensity 
profile (Stokes I). For rapid rotating cool stars the temperature 
information encoded in line profile bumps of the Stokes I pro- 
files is therefore only weakly affected. The dependence of ZDI 
on the temperature inversion however is strong, which means the 
temperature information obtained from the DI must be part of 
the magnetic inversion as the temperature leads to a non-trivial 
scaling of the Stokes V profile amplitude and width. A consecu- 
tive DI-ZDI approach is therefore equivalent to the simultaneous 
inversion as done in this work. 



(30) 



The average sum-of-squares error for a particular line profile can 
then be written as 

(31) 

The mean error E of all individual spectral lines is then given as 

* - k 2> = k 2>M • 

* *=1 * k=\ 



(32) 



On the other hand, if we perform the averaging process before 
minimizing the error function, as we did in our line profile ex- 
tractions, we can write for the mean observed profile pattern 



N k 



(33) 



k=l 



4.2. Multi-line inversion 

For DI as well as for ZDI we could apply the least-squares min- 
imization of Eq. (22) to all of our available and reconstructed 
spectral lines simultaneously, i.e. using a combined long spec- 
tral line vector such that the error function reads 



and the same holds for the mean synthetically calculated profile 



* k=\ 



(34) 



j Np N k N m 

2 2 Z Z (0 "' k ' m - 7 < x w) 2 . 

p=l k=\ m=l 



(29) 



The mean error E for the pre-averaged line profile can then be 
written as 



where x is the model parameter vector, p the rotational phase, 
k the spectral line, m the wavelength or the velocity index re- 
spectively, and N p , Nk, N m their respective total numbers. This 
is the approach usually pursued in DI with TempMap (e.g. Rice 
& Strassmeier 2000; Rice, Strassmeier & Kopf 2011) but fails 
for ZDI since the polarimetric signals are usually buried deep 
within the noise and therefore require a signal extraction which 



E = & 



l k T,(h-Ok) 

* k=l 

This is equivalent to 

j Nt 

N k ^ ek 



(35) 



E = 6 



k=\ ) 



(36) 
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We are assuming again that the noise e k is uncorrelated and of 
zero mean, so that all covariance terms disappear. Then we may 
write 

ly k k=l 

Comparing this with Eq. (32) shows immediately that the mean 
error for the averaged line profile is lower than for the individual 
profiles, 

y k k=l K k=\ 

Most likely, the idealized assumption of completely uncorrelated 
noise is not valid but still, the error of the averaged modeling E 
is not expected to be larger than that of a conventional mini- 
mization of E. Note, that an important aspect here is the correct 
synthetic modeling of the averaged profile. The extension to a 
weighted averaging is straight forward. 

4.3. Differential rotation inversion 

In this section, we introduce our methodology for retrieving a 
surface differential rotation signal, if present. Our technique fol- 
lows closely that introduced by Petit et al. (2002). Given that the 
observations span over several rotational phases and that the stel- 
lar surface obeys a solar-type differential rotation law, the under- 
lying idea is that the closer the modeled differential rotation is to 
the true value the lesser the resulting reconstructed surface image 
is smeared out. In other words, not accounting for differential 
rotation introduces a blurring effect in the resulting DI or ZDI. 
Differential rotation provides not only constraints on the image 
space in terms of the resulting sharpness of the image but also on 
the data space because the differential rotation alters the velocity 
distribution of the individual local line profiles over the visible 
surface and thus produces a modification of the disk-integrated 
flux profile (see, e.g. Gray 2005). The data space approach has 
been extensively used in the past for determining the differen- 
tial rotation of a number of stars (e.g. Reiners & Schmitt 2003; 
Reiners 2006). Compared to the approach chosen by Petit et al. 
(2002) and Donati et al. (2003) there is one major difference we 
apply here. Instead of repeating the inversion many times with 
various differential rotation parameters to sample the x 1 land- 
scape, we directly incorporate the differential rotation law into 
the inversion. All other restrictions given in Sect. 5.1. of Petit 
et al. (2002) also hold for our approach with particular emphasis 
on the required rotational coverage, note that the requirements 
for the phase coverage is less restrictive. 

4.3.1. Assumed differential-rotation law 

The assumption here is that the differential rotation on V410 Tau 
can be described, as for so many cool and late-type stars, with a 
simplified "solar-type" differential rotation law; 

£1(1) = £l eq - dQ sin 2 / , (39) 

where I is the surface latitude, Q(Z) the respective rotation rate at 
a particular latitude, Q. eq the rotation rate at the equator, and <f Q 
the difference in rotation rate between the pole and the equator. 



4.3.2. Implementation as a sheared-image analysis within 
'Map 

The search for the least blurred image consistent with the mini- 
mum of the objective Eq. (22) can also be expressed as an image 
deconvolution problem with a latitudinal smearing operator or 
point-spread function P] at . This concept will be investigated and 
developed in a future paper, instead we concentrate here on de- 
scribing the problem in the same way as in Petit et al. (2002). 
Providing sharpness as an additional constraint in image space 
may indeed provide a means of finding a solution which best 
matches with a differential rotation law and the data error (fit) 
at the same time. This extra constraint has been provided in 
terms of minimizing the information content (maximizing the 
image entropy) of the image by Donati et al. (2000) and Petit 
et al. (2002). However, because maximum entropy is entirely 
invariant under a random permutation of the image data (pixel) 
and does not account for any spatial coherence in the image do- 
main, we believe this constraint is not an appropriate measure to 
quantify distortions due to surface differential rotation. More- 
over differential rotation may, despite the distortions, retain to a 
large extent the morphological coherency of surface structures 
over several rotation periods. Instead of adding an additional 
constraint to the error function, we want to completely rely on 
the achieved quality of the data fit and implement the optimiza- 
tion for the two parameters Q. eq and dQ. within the iMap code. 
This extra flexibility amounts to an additional degree of freedom 
in the time domain to account for the time-dependent rotational 
evolution of the image. The profiles and the corresponding part 
of the visible stellar disk are no longer only determined just by 
the spectral information of the respective phase alone but also 
by the time information (Petit et al. 2002). This minimization 
process can conveniently be incorporated into our Landweber it- 
eration in Eq. (28) to yield an alternating minimization scheme. 
If we define the vector £2 which comprises our two differential 
rotation parameter £l eq and dQ. the process can be described as 
follows: After proceeding along the negative gradient direction 
in image space (temperature and/or magnetic field) the process 
switches to the omega space to perform a minimization of the 
squared error, Eq. (22) with respect to the differential parameter 
vector il. This provides a new current estimate of the vector £1. 
which then enters again in the next iteration cycle of Eq. (28) 
This process is repeated until the stopping condition is reached 
(see Sect. 4). We can formally introduce an image-shear opera- 
tor P d which differentially rotates a current estimate x k (il) of the 
image (temperature or magnetic field) to write the process as 

Q = argmin ¥ d [x k (£l) + w k V(x k (£l)) (O - I(x*(tt)))] , (40) 

The minimization of Eq. (40) is performed by a gradient descent 
method similar to the one used for the Landweber iteration. 

4.4. Fixed stellar parameters 

The adopted stellar parameters are the same as those used and 
partly derived in Paper I. The effective temperature is 4,400 K, 
the logarithmic gravity is 4.0, the projected rotational velocity 
is 77.7 kms -1 , and the inclination of the stellar rotation axis is 
62° as determined from DI in Paper I. The micro- and macro- 
turbulence is set to 2.0 kms -1 and 1.5 kms -1 , respectively, so- 
lar chemical abundances are assumed for all elements. As also 
described in detail in Paper I, an average rotational period of 
1.871948d was adopted. Rotational phase is also kept the same 
as in Paper I. 
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Fig. 5. Stacked plot of the synthetic fits (thick dashed lines) to the observed Stokes profiles (thin solid lines), on the left side Stokes I and on the 
right Stokes V. All spectral line profiles are normalized by the continuum intensity. 



4.5. iMap setup 

For the inversion with iMap, we use all 12 available Stokes I 
and Stokes V SVD profiles from Fig. 3 and Fig. 4. We utilize 
the multi-line approach described in Sect. 4.2 which requires to 
synthesize all 929 spectral line profiles of the observation ma- 
trix to finally obtain the mean synthesized profile. The enor- 
mous amount of computation is performed by our ANN ap- 
proach (Carroll, Kopf & Strassmeier 2008). As mentioned in 
Sect. 4. 1 all line parameters were taken from the VALD line list 
and model atmospheres are produced by ATLAS-9 (Castelli & 
Kurucz 2004). The stellar surface is partitioned into 6° x 6° seg- 



ments, whereas the polarized radiative transfer is calculated on 
a subpartition of 3° x 3°. The initialization temperature was ho- 
mogeneously set to 4,400 K. The initial magnetic field is homo- 
geneously set to zero for all three magnetic components (radial, 
azimuthal and meridional). For v sin i, stellar inclination, gravity, 
and abundances, the values given in Sect. 4.4 are used. 

The stopping rule according to the discrepancy principle for 
the Landweber iteration is set according to the standard errors 
of the reconstructed line profiles. As described in Sect.4.1 this 
means that we use for the magnetic inversion the largest (maxi- 
mum) standard error of the reconstructed Stokes V profiles and 
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for the temperature inversion the largest standard error of the re- 
constructed Stokes I profiles. The individual values are shown in 
Fig. 3 and Fig. 4. 

5. Application to V410 Tau 

5.1. Temperature distribution in 2008/09 

The synthetic fits relative to the observed Stokes I profiles are 
shown in a stacked plot on the left side of Fig. 5. The final tem- 
perature map in orthographic projection is shown in Fig. 6. It is 
dominated by a structured, cool, polar cap with a temperature as 
cool as «3,500 K, i.e. »900 K below the effective temperature. 
Several appendages reach down to lower latitudes of around 30° 
with decreasing temperature contrast the further they are away 
from the cap. The total spot filling factor is 11.9 %. From pre- 
vious investigations, we knew that V410 Tau shows a polar spot 
and a number of smaller cool spots together with localized hot 
features (see the review in Paper I). However, the simultaneous 
reconstruction of the surface temperature and the magnetic-field 
distribution provides unprecedented accurate surface details. 

An interesting feature of the temperature map in Fig. 6 is the 
pair of isolated spots consisting of a cool spot at phase <p « 0.7 
and of a hot spot at phase <p « 0.6 at intermediate latitudes. The 
hot spot's temperature is «4,800 K, i.e. »400 K above the effec- 
tive temperature while the cool spot's temperature is «4,000 K, 
i.e. »400 K below the effective temperature. A similar pair of 
spots, although larger in size and with much smaller contrast, is 
located on the equator at phase (p ~ 0.25. This pair's temperature 
difference is a400 K in total. Yet a third such pair of just mod- 
erate warm and cool spots may exist around <p « 0.0 at latitude 
~ 30°. 

All of above features were independently reconstructed with 
TempMap and iMap in Paper I; some with slightly different con- 
trast but the overall similarity is very good. 

5.2. Magnetic field distribution in 2008/09 

The magnetic field reconstruction found by iMap provides a so- 
lution compatible with the observed data and the noise level. The 
synthetic fits for the individual Stokes V profiles are shown in a 
stacked plot in Fig. 5. The ZDI map is shown as an orthographic 
plot in Fig. 7. The distribution of the individual magnetic com- 
ponents (radial, azimuthal and meridional) can be best seen in 
the pseudo Mercator-projection plot in Fig. 9. The orientations 
of the field within the Mercator plots is defined such that a pos- 
itive direction in the radial field points outward while a negative 
polarity points inward, the meridional component is positive if it 
points toward the upper north-pole and negative in the direction 
of the south-pole. For the azimuthal component a positive field 
is defined in westward direction and negative for an eastward 
direction. 

The magnetic topology on V410 Tau is dominated by a bipo- 
lar structure around the visible pole with peak field strengths of 
up to ±1.9 kG. Figure 8 emphasizes the relatively good corre- 
lation of the temperature and magnetic field and shows that the 
strongest fields are located in the darkest regions around the pole. 
The two polarities appear separated by a sharp, pole-crossing, 
neutral sheet. In total the field forms an S-shaped structure cen- 
tered at the rotation pole. The topology is predominantly radial 
but 1 kG meridional components coexist very close to the pole. 
An azimuthal component, on the other hand, is restricted to a re- 
gion at «60° -latitude and around 270° longitude, i.e. 180° away 
from the strongest radial component, close in distance but on the 



opposite hemisphere and always of positive polarity. Its peak 
field strength reaches 1 kG around 270° longitude while the rest 
of it mostly remains near 500 GOne can express the field topol- 
ogy also in terms of the poloidal and toroidal components which 
reveals in our case that the majority, 73 %, of the surface field of 
V410 Tau is given in the form of poloidal fields and just 27 % 
being in form of a toroidal component. 

The dominant of the two isolated warm-cool spot pairs 
around <p « 0.6 and <p ~ 0.7 has no recognizable distinct 
magnetic-field structure. Both, the cool and the warm spot, ap- 
pear with the same general mix of field components as their im- 
mediate surrounding and at least this particular spot pair is re- 
constructed with the same (positive) polarity. On the contrary, 
the second such spot pairs (around (p as 0.25; see Sect. 5.1 above) 
indicates a mixed polarity, negative for the warm feature and 
positive for the cool feature. 

It is interesting to note that the strong Ha flare reported in 
Paper I (located at <p ~ 0.18) coincides with the best visibility 
of the sharp neutral sheet of the bi-polar radial field. It is at 
least tempting to state that, if truly coincident, it would very well 
agree with solar active-region physics where flares are usually 
related to the current sheets separating opposite polarities (e.g. 
Schrijver et al. 2005). 

An error analysis of the here presented results is given in the 
appendices A, B, and C. 

5.3. Differential rotation of V410 Tau 

As described in Sect. 4.3 a solar like differential rotation param- 
eterized by Q. eq and dO. is part of the inversion process. In the 
current version of iMap the sheared image analysis is applied 
to the temperature structures only. For the initialization, Q. eq is 
derived from the rotational period of 1.871948 d (see Paper I). 
Only a very weak solar-like differential rotation is present, if at 
all. The formal value for the angular velocity at the equator is 
£l eq = 3.356±0.005, well within the errors of the photometrically 
determined rotation period, and for the differential rotation rate 
between the equator and the pole we have c/Q = 0.007 ± 0.009. 
We can not safely regard this as a detection of a differential rota- 
tion on V410 Tau but conclude that it must be very small, likely 
smaller than a factor 30 compared to the Sun. 

5.4. The magnetic nature of V410 Tau 

The probable absence of a radiative core and the non-detection of 
differential rotation in this work and other investigations (Strass- 
meier, Welty & Rice 1994; Skelly et al. 2010) leaves not much 
room for speculation about the nature of the dynamo operating in 
the interior of V410 Tau - at least in the framework of classical 
mean field dynamos. As has been shown by Kiiker & Rudiger 
(1999) a likely explanation for the generation of the magnetic 
field in weak-lined T Tauri stars can be given in terms of a a 2 - 
dynamo. 

From the observational point of view we can characterize 
and quantify the surface topology with various degrees of de- 
tail but providing a quantitative link between observation and 
theory is anything but straightforward. In the moment we feel 
that dynamo theory does not give enough reliable information to 
provide such a quantitative link between the surface appearance 
of magnetic fields and their internal generation process. But this 
shall not keep us from giving qualitative arguments in favor of an 
a 2 -dynamo here. Following Kiiker & Rudiger (1999) for an a 2 - 
dynamo the field topology exhibit a distinct non-axisymmetric 
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Fig. 6. Orthographic maps of the temperature distribution of V410 Tau at four different rotational phases <p. 




Fig. 7. Orthographic maps of the magnetic field distribution of V410 Tau at four different rotational phases (/>. The magnetic field is color coded 
and the field lines are proportional to the field strength (blue negative and red positive polarity). 




Fig. 8. Orthographic maps of the temperature overplotted by the magnetic field lines at four different rotational phases <p. For a better visibility 
the positive fields are color coded in yellow while the negative fields are again in blue. As above the field lines are proportional to the strength of 
the field. 



geometry with SI- or A 1 -type modes (Riidiger & Hollerbach 

2004) . The pronounced bipolar magnetic field at the pole (Fig. 
7) can easily be identified with such a non-axisymmetric ge- 
ometry. However the fact that the non-axisymmetric mode of 
such an a 1 -dynamo solution can only rotate with the star itself 
(except for a slight drift) may not really fit to recent observa- 
tions that V410 Tau exhibit a periodic change of its activity in 
the range between 4.8 to 5.5 years (Stelzer et al. 2003; Savanov 
2012). But a possible explanation can still be given in terms of 
an a 2 f2-dynamo (Elstner & Korhonen 2005), where the non- 
axisymmetric solution is modulated by an oscillatory axisym- 
metric mode. The excitation of such an ff 2 £2-dynamo can already 
occur at a very low differential rotation rate (Elstner & Korhonen 

2005) . To qualitatively compare this dynamo model with our re- 
constructions we have taken a snapshot of the non-axisymmetric 



mode of these simulations at a time when both magnetic polar- 
ities are of approximately the same strength. The surface field 
strength was scaled to be comparable to that of the reconstructed 
surface field of V410 Tau and the symmetric field (relative to the 
equator) at the south pole has been damped. The similarity be- 
tween our reconstruction and the o^fi-dynamo simulation shown 
in Fig. 10 is quite remarkable and may in fact strengthen the hy- 
pothesis that the field of weak-lined T Tauri stars is generated by 
an a 2 -type dynamo or in cases where a weak differential rotation 
has already developed by an or 2 £2-dynamo. 

An interesting model which accounts for the evolution of 
pre-main sequence dynamos were presented by Kitchatinov et 
al. (2001) and may also be applicable to our findings. Accord- 
ing to their simulations the dynamo during the early phases of 
the pre-main sequence stage, when the core development is still 
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Fig. 9. Pseudo-Mercator projection of the three magnetic-field components; left radial, middle meridional, and right azimuthal. Note, that the 
lower x-axis is the longitudinal angle, while the upper x-axis gives the rotational phase. The dotted line indicates the limit for the visibility due to 
the stellar inclination. 



in progress and the rotational velocity is large, gives rise to a 
strong non-axisymmetric mode. In the course of the evolution 
when differential rotation sets in the star develops an increas- 
ing axisymmetric dynamo mode which can lead to an oscillating 
behavior (Riidiger & Elstner 1994) before the field is eventually 
completely dominated by the axisymmetric mode at a later stage. 

How does V410 Tau compares with other T Tauri stars that 
have been subject to magnetic investigations. Zeeman-Doppler- 
Imaging has been applied to a small number of classical T Tauri 
stars and have been compared by Hussain (2012) to find a re- 
lation between mass, period, magnetic field strength and field 
complexity. Although no clear relation could be found from the 
limited sample, there seems to be a tendency that fully convec- 
tive T Tauri stars may harbor simpler surface fields (i.e. mainly 
a dipole pole component) while T Tauri stars which have devel- 
oped a radiative core have the tendency to more complex fields 
Hussain (2012). It is not clear if V410 Tau as a weak-lined T 
Tauri fit into this picture. Given its evolutionary stage it may 
be reasonable to assume, in the light of the above trend, that 
V410 Tau exhibit a more complex surface field than just a dipole 
(spherical harmonic numbers I = 1 and m = 0). This is in fact 
true, but our magnetic reconstruction shows just a little more 
additional complexity with a clear dominance of the SI mode 
(I = 1 , m = 1) and the Al mode (/ = 2, m = 1). The apparent 
diversity of the magnetic field of T Tauri stars certainly requires 
a greater sample in order to provide a more general picture of 
the magnetic evolution and the underlying dynamo for these pre- 
main sequence stars. 

At the end we also want to comment on the apparent dif- 
ferences between the reconstructed surface fields of Skelly et 
al. (2010) and our magnetic and temperature maps. In a re- 
cent magnetic investigation of V410 Tau Skelly et al. (2010) 
found that 50 % of the total field has a toroidal component 
and the other 50 % are poloidal and that the majority of the 
field shows strongly inclined large-scale azimuthal fields which 
are distributed over one entire hemisphere, whereas our recon- 
structed topology has only 27 % toroidal fields and 73 % poloidal 
fields with a strong radial and bipolar component around the 
polar region. Moreover, in contrast to our reconstructed maps 
where the strong magnetic fields are well correlated with the 
temperature, see Fig. 8, the fields of Skelly et al. (2010) seem 
to show almost no correlation with their brightness maps (i.e. a 
proxy of the temperature distribution). The fact that these data 
were observed at almost the same time as ours (January 2 to Jan- 
uary 17, 2009) at the Telescope Bernard Lyotl (TBL) with the 
NARVAL spectropolarimeter at the Pic di Midi, is certainly a 



reason for a closer inspection and comparison of the two results. 
Without having the possibility for doing such a detailed compar- 
ison of the data as well as of the inversion codes we can only 
speculate about the reasons for the disparity of the two recon- 
structions. But we want to emphasize that in contrast to Skelly 
et al. (2010) our approach makes no assumptions on the surface 
topology in terms of a spherical harmonic decomposition or po- 
tential field structure. It is furthermore fully based on polarized 
radiative transfer and line profile modeling and it pursues a strat- 
egy which simultaneously invert the temperature and magnetic 
field. This give us enough confidence to believe in the validity 
of our results. 



6. Summary and conclusions 

Despite the relative high levels of magnetic activity of many 
rapidly rotating late-type stars the typical spectropolarimetric 
signals that these stars produce are mostly below the noise level 
of a typical spectropolarimetric observation and require there- 
fore specialized techniques to extract and reconstruct the line 
profiles. To address this problem we have developed a new SVD 
based spectral-line extraction technique. For Stokes V this ap- 
proach is a full reconstruction method which allows us in this 
work to reach an improvement in the signal-to-noise ratio by a 
factor of 40. For Stokes I where the actual signals are orders of 
magnitudes larger then for Stokes V the SVD method reduces to 
a denoising technique but still allows to increase the signal-to- 
noise ratio by a factor of 5. 

The interpretation of polarimetric data in terms of ZDI al- 
ways require the introduction of additional information. This 
is true for the underlying model used in the inverse problem as 
well as for the regularization of the problem itself. To pose a 
minimum number of surface constraints on the solution we have 
chosen a new regularization approach for DI and ZDI in terms 
of an iterative regularization of the problem. The new variant 
implemented in the inverse module of iMap shows a good and 
stable behavior and converges within 400 iterations. 

A common problem for ZDI is that the real strong magnetic 
fields which are associated with cool surface regions produce 
only a fraction of the photon flux compared to the unaffected 
quiet or even hot surface regions. To prevent our code from 
these ill-defined flux-weighting and the mutual effects between 
the temperature and magnetic fields we have pursued a strategy 
that performs a simultaneous DI and ZDI to retrieve both the 
temperature and the magnetic field distribution at the same time. 



Article number, page 13 of 18 




Fig. 10. Intriguing structural similarity between the result of our magnetic surface reconstruction (left) and an a 2 fl-dynamo simulation which is 
scaled to a surface field strength of 1.9 kG (right). 



The reconstruction of the temperature and magnetic field of 
V410 Tau reveals that the majority of the strong fields of 2 kG 
are located in the cool spots, in particular within the large polar 
spot. The reconstruction clearly show that both polarities coexist 
within the large polar spot and that the entire polar-field topology 
appears to be dominated by a twisted bi -polar structure separated 
by a relatively sharp neutral line. It is reminiscent of an over- 
dimensioned solar active region. Even the time of a strong flare 
coincides with the best visibility of this neutral line, just like for 
solar active regions. However, due to its shear size of approx- 
imately 11% of a hemisphere, the bi-polar feature appears S- 
shaped as if it had been dragged around the pole due to the rapid 
stellar rotation. The absence of a detectable differential rotation 
and the pronounced non-axisymmetry of the field may suggest 
an a 2 -dynamo operating in the interior of V410 Tau. This is also 
supported by the intriguing similarity between our reconstructed 
surface fields and the dynamo simulation of Elstner & Korhonen 
(2005), and may support the theory that the underlying mecha- 
nism responsible for the magnetic field generation in weak-line 
T Tauri stars is an cr 2 -type dynamo. 
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Appendix A: Error Analysis 

In this section we want to investigate the significance of our re- 
construction results, in particular we want to quantify how sen- 
sitive the solution of the inverse problem is with respect to the 
initial conditions. In our DI and ZDI setup we deal with four 
parameter spaces (temperature, radial magnetic field, azimuthal 
magnetic field and meridional magnetic field) each of which with 
a dimension equal to the number of surface segments. The inver- 
sion algorithm has to navigate through the combined parameter 
spaces to find a solution that is compatible with the data. In order 
to study the stability of the inversion relative to the initial start- 
ing conditions (different positions in the parameter spaces) we 
use a simulation that runs the inversion with the original data set 
but from randomly chosen starting positions. Each parameter 
space is independently initialized by choosing a random value 
for each surface segment. Though we may create a consider- 
able spread among the individual surface segments we can only 
consider a small fraction of the overall parameter space. Any ex- 
haustive study to sample the parameter space is way out of reach 
for this high-dimensional setting and it is not the purpose here, 
instead we want to gain a quantitative measure how the result 
varies when the initial conditions are changed. 

A random generator on the basis of a normal distribution pro- 
vides the values for each parameter space and surface segment. 
The underlying normal distribution has a mean value which is 
set to the effective temperature for the temperature component 
and to zero Gauss for all three magnetic components. The stan- 
dard deviation is 250 K for the temperature and 250 G for all the 
magnetic components to provide a large spread among the indi- 
vidual surface segment values. The simulation (inversion with 
different initializations for all parameters and segments) is run 
for 100 times on the original data set. For the so obtained 100 
DI and ZDI maps we calculated the root-mean-square (RMS) 
error to finally get RMS error maps for each parameter. The in- 
dividual inversions are all run to the same accuracy level given 
by the noise of the reconstructed data. Since the derivatives with 
respect to the magnetic parameters as well as with respect to the 
temperature become to weak below a latitude of -30 degrees to 



provide any substantial changes during the inversion we fix these 
segment values to zero for the magnetic maps and to the effective 
temperature for the temperature maps. 

The RMS error map for the temperature is shown in the or- 
thographic plot Fig. A.l. It shows a peak value of 80 K at po- 
sitions associated with the coolest regions (polar spot). Other 
regions hardly exceed a RMS error of 50 K. The RMS error 
map for the three magnetic components is shown in Fig. A. 2, 
Fig. A. 3, and Fig. A.4. The error map of the radial magnetic 
field has a peak value of 125 G, the average RMS error value is 
64 G. The RMS error values for the meridional and azimuthal 
field exhibit peak values of 101 G and 106 G respectively and 
both have a mean RMS error value around 60 G. What can be 
seen from these maps is that the error values are correlated with 
the absolute field strength as well as with the temperature. This 
emphasizes once more the influence of the temperature on the 
magnetic field determination. For the given photospheric condi- 
tions of V410 Tau a temperature change of just 80 K in a spot 
like region causes a change in the amplitude of the Stokes V 
signal of 4%. In the strong field regime of the polar spot this 
difference in amplitude is equivalent to a magnetic field of 60 G 
! Given this temperature dependence, the noise level in the data 
as well as the fact the a random small scale magnetic field on the 
surface produces local Stokes V signals that effectively cancel 
each other out, the obtained uncertainties from the error simu- 
lation are surprisingly low and show that the inversion always 
settle in the proximity of the same (local) minimum and that the 
solution is robust against perturbations of the initial conditions. 



Appendix B: Test Inversion with a dynamo model 

In this section we want to examine if the inversion is able to 
retrieve a surface field topology as complicated as that of a sim- 
ulated ff 2 n-type dynamo with the given data set. The question 
that we want to investigate here goes beyond the testing of an 
inversion code with simplified surface topologies, what we want 
to do here is to seek a more direct answer to our problem, can 
we trust the hypothesis of our work that the reconstructed sur- 
face topology of V410 Tau is in fact similar to that of a given 
a 2 0-type dynamo simulation. 

For that reason we synthesized from the a 2 £2-dynamo simu- 
lation of Elstner & Korhonen (2005), 12 phase resolved spectra 
at phases which correspond to our observations. The dynamo 
model was scaled to have a peak field strength of 1920 G. The 
stellar parameter are those from V410 Tau given in Sect. 4.4. 
Noise is added according to the Bootstrap analysis in Sect. 3.2. 
The inversion with iMap is run with the same setup as the origi- 
nal inversion for V410 Tau, which also means that the inversion 
is stopped at the same noise level. Only Stokes V profiles are 
used for this magnetic inversion. The result is illustrated in the 
orthographic plots Fig. B.2 and Fig. B.3 as well as in the profile 
plot Fig. B. 1 . Given the number of available phases and the noise 
level in the data, the original magnetic field structure is remark- 
ably good reproduced. Although the maximum field strength is 
not quite reached in the reconstruction shown in Fig. B.3 one 
can see that the intertwined nature of the field structure and the 
rapid variation of the field lines around the polar region is well 
reproduced. The quality of the fit of the inversion relative to the 
original synthesized Stokes V profiles is shown in Fig. B.l. This 
demonstrates that with the given data set the inversion code can 
in fact reconstruct a surface field with a topology as complex as 
a a 2 Q-type dynamo. 
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Fig. B.l. Stacked plot of the synthetic fits (dotted points) and the 
synthetic Stokes V profiles (grey solid lines) obtained from the dynamo 
model. 



Appendix C: The influence of line blends 

Here we want to investigate how much information about the line 
blends may get lost within the rank estimation of the truncated 
SVD procedure. Note that in our analysis we compare the SVD 
reconstructed line profile with the synthetic (weighted) mean of 



all contributing lines (929 lines for Stokes V and 56 for Stokes 
I) where all significant line blends are accounted for. One might 
suspect that if information about the line blends is leaked into 
the noise subspace it may lead to a mismatch between the syn- 
thesized line profiles and the reconstructed profile which may 
cause problems in the subsequent inversion. 

Let us concentrate on the Stokes V profiles here and have a 
look on how line blends influence the resulting profile from the 
two-stage SVD reconstruction. Line blends affect the individual 
lines used in the observation matrix in different ways, they occur 
at different positions within the line profiles, and are of different 
strength and number. So it will be of interest to see how much 
of the blends are finally recognized as systematic and correlated 
effects within the SVD reconstruction. For that purpose we have 
modeled a simple monopole star with a homogeneous magnetic 
field of 10 G and a small rotational velocity of 5 km s _1 . All 
other stellar parameters are assumed to be the same as that for 
V410 Tau. We synthesize the line profiles in the velocity do- 
main and use a value of + 150 km s _1 , around the respective line 
center of each contributing line. For such a small rotational ve- 
locity the Stokes V signals are therefore much thinner than the 
used velocity range. This wide span of the velocity domain will 
give us additional information about the systematics introduced 
by the line blends relative to the base level (i.e. zero polarized 
continuum). The 929 lines are calculated without any noise. Af- 
ter calculating the line profiles and their respective line blends 
we applied the SVD analysis and used all available eigenprofiles 
for the reconstruction. The result of the two stage process is an 
average Stokes V profile that is the same as a regular average 
since we have used the entire set of eigenprofiles. In Fig. C. 1 on 
the right side we see the result of the SVD reconstruction. It is 
clearly seen that the Stokes V profile has its regular shape in spite 
of the many contributing line blends. The wiggling is the result 
of the line blends. Note again no noise is present in that synthetic 
test case. This wiggling of the base line (the continuum) shows 
no dramatic variation within the velocity range. We may quan- 
tify this wiggling in terms of a noise measure, and use the me- 
dian absolute deviation (MAD) for that purpose which is defined 
as, MAD(V) = mediani\V\ - median(Y)\, \V„ - median(Y)\), 
where median(V) is the median of the vector components of V. 
The MAD give us a value of 1.56 x 10 -5 . Given the strength of 
the reconstructed Stokes V signal, which has an amplitude com- 
parable to that observed for V410 Tau, we can see that the blends 
introduce an effect that is below the noise level of the real recon- 
structed observations. How much information is lost if we would 
just use the eigenprofile corresponding to the largest eigenvalue? 
In that case we would pretend that our signal subspace is of di- 
mension one and all significant signal information is comprised 
within the first eigenprofile. Performing the two-stage SVD re- 
construction under this assumption results in a Stokes V profile 
shown on the right side in Fig. C.l. The difference between the 
reconstruction with the full set of eigenvectors and that using 
only a single eigenprofile is very small. We calculated the RMS 
between the two reconstructed profiles which gives a value of 
2.22 x 10 -5 . The contribution of the line blends relative to the 
Stokes signal for the set of 929 spectral lines used in this work 
is apparently very small, i.e. ten times smaller than the noise 
level deduced in Sect. 3.2. The remaining systematic effects 
introduced by the blending are essentially captured by the first 
eigenprofile which is also the signal carrying component. This 
demonstrates that the influence of line blends for the set of spec- 
tral lines used in this work is small and that the leakage of infor- 
mation carried by eigenprofiles belonging to smaller eigenvalues 
(i.e. noise space) is negligible. 
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Fig. A.l. Temperature error map in orthographic projection at four different rotational phases. The maximum RMS error value is 80 K and the 
mean RMS error value is 50 K. Regions below a latitude of -30° are set to a constant value and do not participate in the error estimation. 
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Fig. A.2. Radial magnetic field error map in orthographic projection at four rotational phases, the maximum RMS error value is 125 G and the 
mean RMS error value is 64 G. 
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Fig. A.3. Meridional magnetic field error in orthographic projection, the maximum RMS error value here is 101 G and the mean RMS error 
value is 5 1 G. 
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Fig. A.4. Azimuthal magnetic field error in orthographic projection with a maximum RMS error value of 106 G and a mean value of 52 G. 



Article number, page 17 of 18 



<p= 0.00 cp= 0.25 <p= 0.50 <p= 0.75 




Fig. B.2. Orthographic maps of the original dynamo simulation which serves as a synthetic model for the test inversion, The model is shown at 
four different rotational phases (/>. 
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Fig. B.3. Reconstructed surface field from 12 phase modulated Stokes V profiles. The number and order of the phases as well as the introduced 
noise correspond to the real observations of V410 Tau. Given the observational constraints imposed on this test inversion the dynamo model could 
be very well reproduced. 
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Fig. C.l. Reconstructed SVD profile originating from all 929 spectral lines (left). The reconstruction is performed with all eigenprofiles. The 
"quasi-noise" level introduced by the line blends has a median absolute deviation (MAD) of 1.56 x 1CT 3 . On the right, reconstructed SVD profile 
with only one eigenprofile corresponding to the largest eigenvalue. The difference between the reconstruction using the full set of eigenprofiles is 
hardly visible and the RMS between them is 2.22 x 1CT 5 . 
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